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The hydrodynamical description of the "Little Bang" in heavy ion collisions is surprisingly suc- 
cessful, mostly due to the very small viscosity of the Quark-Gluon plasma. In this paper we sys- 
tematically study the propagation of small perturbations, also treated hydrodynamically. We start 
with a number of known techniques allowing for the analytic calculation of the propagation of small 
perturbations on top of the expanding fireball. The simplest approximation is the "geometric acous- 
tics" , which substitutes the wave equation by mechanical equations for the propagating "phonons" . 
Next we turn to the case in which variables can be separated, where one can obtain not only the 
eikonal phases but also the amplitudes of the perturbation. Finally, we focus on the so called Gubser 
flow, a particular conformal analytic solution for the fireball expansion, on top of which one can 
derive closed equations for small perturbations. Perfect hydrodynamics allows all variables to be 
separated and all equations to be solved in terms of known special functions. We can thus col- 
lect the analytical expression for all the harmonics and reconstruct the complete Green function of 
the problem. In the viscous case the equations still allow for variable separation, but one of the 
equations has to be solved numerically. Summing all the harmonics we show real-time perturbation 
evolution, observing the viscosity-induced changes in the spectra and the correlation functions. The 
calculated angular shape of the correlation function is remarkably similar to the shape emerging 
from the experimental data, for sufficiently large viscosity. We predict a minimum at m ~ 7 and 
maximum at m ~ 9 harmonics, which also have some experimental evidence for it. We conclude 
that local "hot spots" in the initial state are the only visible origin of the observed correlations. 

PACS numbers: 



I. INTRODUCTION 

Since it is the third paper of the series devoted to the 
propagation of perturbations on top of the "Little Bang" , 
it does not need a detailed introduction. Let us only 
briefly point out the main physics of the phenomena in 
question, and then mention where the reader can find 
important earlier works on the subject. 

Initial state perturbations of an "average fireball", 
which occur on an event-by-event basis, lead to diver- 
gent sound waves, similar to the circles from a stone 
thrown into a pond. The sound velocity is ^ 1/2 and 
the time till freezeout Tpo ~ 2i? (where R is the nuclear 
size, about 6 fni for Au nuclei used in the experiment), 
thus the "sound horizon" (the maximal radius of the cir- 
cles) reaches Hg ~ R. In terms of the angular variables 
we use, it means a response at relatively large angles, 
0{zLl radian), from the perturbation. The strong radial 
explosion of the fireball dramatically enhances the con- 
trast, making small deviations of the freezeout surface 
easily observable experimentally, provided the transverse 
momenta of the particles are tuned into the appropriate 
range. The shape of the hydro response to an initial point 
perturbation (the Green function) is quite non-trivial, 
and we show that for appropriate values of the viscos- 
ity it reproduces the shapes of the two-point correlation 
functions observed experimentally surprisingly well. We 
will conclude with the "minimal" and "maximally coher- 
ent" scenarios of the collisions: for experimental selection 
between those one needs to measure certain three-point 



correlations functions, as was discussed in detail in our 
previous paper [T]. 

Many issues we discuss, such as the power spectrum of 
higher harmonics of perturbations, are analogous to the 
events in Cosmology during the last decade. We mean 
in particular the observations of the sound horizon scale, 
both in the cosmic microwave background (CMB) radi- 
ation (see e.g. [2] and the earlier work cited in it) and 
in the distribution of galaxies [3] . Discussing similarities 
and differences between the Little and Big Bangs will be a 
recurring theme of this paper. Let us just comment that 
while these observations did turn Cosmology into a much 
more quantitative science, hopefully their "Little Bang" 
analogues will also help us to fix the global parameters 
of nuclear collisions and the QGP much better. 

Outlining the paper's context, we now go into a bit 
more detail over the brief history of the "second act of 
hydro" . Sound propagation on top of the expanding fire- 
ball was first considered by Casalderrey-Solana and one 
of us (ES) in [3]. The fireball expansion was modelled 
by a Universe expansion using the Friedmann-Lemetre- 
Robertson- Walker metrics, and the specific phenomena 
discussed in it was the effect of the variable speed of 
sound (due to the QCD phase transition) on sound prop- 
agation. Its main result was the appearance of backward- 
moving or convergent spherical/conical waves, together 
with the usual divergent ones. It is worth noting that the 
hadronic era has a near-constant speed of sound cs ~ 0.4 
as noted in |5] and established later for the chemically 
non-equilibrated version of hadronic matter in p| . The 
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so called mixed phase era is the only one in which Cg 
varies. 

A qualitative picture of the "sound circles" resulting 
from point-like initial-state perturbations, and reaching 
by freezeout the so called "sound horizon" radius, were 
introduced in the first paper of this series [?]• 

(It is amusing that Gurzadyan and Penrose [5] not only 
came claim that the WMAP data provide some evidences 
for circles, or even co-central circles, in the CMB tem- 
perature variations. The ones they found, however, have 
sizes few times larger than the sound horizon scale. So, 
if the claim is statistically sound, those must be some 
pre-Big-Bang events.) 

Unlike the Big Bang, for which one reads the tempera- 
ture perturbations from the sky, the observable traces of 
the sound circles in the Little Bang are not so direct. The 
temperature and velocity perturbations both contribute 
to the particle spectra at the freezeout, and the picture 
is strongly affected by strong radial flow and the exis- 
tence of the fireball's boundaries. The contribution of all 
of this to the spectra predicted in was the "double- 
horn" shape of the angular distribution, with two max- 
ima identified with the latest crossings of the sound circle 
with the fireball boundaries. The "circle" phenomenon 
has also been found by the Brazilian group, in their (zero 
viscosity) numerical studies of "event-by-event hydrody- 
namics" |9] . This group however went further and calcu- 
lated the two-body correlators, finding their characteris- 
tic three-maxima structure. The details of such structure 
in our (viscous) solution will be compared to the experi- 
mental data at the end of this paper. 

A general setting of the problem, including the iden- 
tification of the two basic scales of the problem, the so 
called "sound horizon" and "viscous horizon" , was made 
in the second paper of the series [T], in which we also 
studied in detail the perturbations using the geometric 
Glauber model. Similar ideas have also been proposed 
by Mocsy and Sorensen in [101 EI] • 

The impetus for experimental studies of perturbation- 
related effects was provided by the paper by Alver and 
Roland [T^. They have pointed out that the two par- 
ticle correlation data contains large third angular har- 
monics, and attributed it to the "triangular" shapes of 
some events. That prompted many studies of the ini- 
tial perturbations in the Glauber model, in which the 
fluctuations are due to the random positions and the in- 
teraction probability of the colliding nucleons inside the 
nuclei. It has been found that the data can indeed be 
explained, by Glauber estimates of the initial perturba- 
tions < £3(171)^ > times the "hydro response" at freezeout 
(fo) V3(/o)/e3(m). 

In general, there are two different views on the na- 
ture of the perturbations. A priori, the structure of the 
initial state perturbations can either (i) be just Gaus- 
sian noise or (ii) contain important correlations between 
the harmonics. In the former case, the "minimal Gaus- 
sian scenario" of the initial state, the set of input pa- 
rameters < em(inY > has all the information one may 



possibly need, and all that needs to be done is the hy- 
drodynamical calculation of the "linear response" ratios 
Vn{Pt, fo)/€n{in) between the initial perturbations of the 
fireball shape and the final fiow, for each of the har- 
monics. The other school was pioneered by the above- 
mentioned Brazilian group, that started to do "event- 
by-event hydrodynamics" for many (hundred thousands) 
initial conditions, provided by certain event generator. 
Clearly this only makes sense if one hopes to reproduce 
certain non-trivial correlations contained in experiment 
in a statistically significant way. 

Our study in the previous work [T] , based on Glauber 
theory, had indeed found non-trivial phase correlations 
between all odd harmonics m = 1,3,5.... We have as- 
cribed those to the so called "hot/cold spots" in the ini- 
tial matter distributions, which can appear at any angle 
and are mutually uncorrelated. We also pointed out the 
role of the higher correlators and the "resonance condi- 
tion" between three (or more) harmonics in order to mea- 
sure the relative phases. Similar studies have also been 
done elsewhere, see e.g focused on the resonance be- 
tween the first and third harmonic with the second (the 
reaction plane) and the triangular flow. 

Our main aim in this work is to derive the magnitude of 
all harmonics of the flow in the same setting. Only then 
one can study their coherent sum, the Green function etc. 
This goal is achieved (semi) analytically, with separation 
of variables and full inclusion of viscosity effects. New 
important phenomenon - existence of acoustic dips and 
peaks in the power spectrum - is suggested, calculated 
and correlated with experimental data. 

However these breakthroughs came with a prize: we 
consider (i) only the central collisions, (ii) only conformal 
EoS of matter and (iii) only small perturbations. The 
reader should be aware of the fact that our results should 
capture the qualitative behaviour rather than produce 
accurate numbers, directly related to the experimental 
data. Corrections to non-conformity and non-linearity as 
well as not-too-large non-centrality can be also studied, 
but those will be done elsewhere. 

Let us only comment here on the issue of non-linearity. 
If (as we believe) all harmonics add up coherently, the 
perturbations are generally not small, 0(1), at the initial 
time. However, as the perturbation expands and becomes 
a large sound circle (with the radius up to the "sound 
horizon" size comparable to that of the fireball itself) , it 
quickly becomes small. This is especially true for higher 
harmonics (to which this paper is mostly devoted) as they 
are additionally suppressed by viscosity. 

Clearly, early time evolution of perturbations is non- 
linear, and releated effects are not captured by our ap- 
proach. One practically important issue here is the speed 
of the waves, which affects the size of the sound circle 
at freezeout, which subsequently determines positions of 
the peaks in the correlation function and the power spec- 
trum. Finite amplitude waves are known to travel faster 
than sound. We had investigated this correction and will 
include it in our subsequent paper. Its effect is rather 
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modest: for example a factor 2 matter compression leads 
to only 15% increase in speed. Note, that realistic EoS 
leads to the speed of sound at late stages of the collision 
to be about 20% lower than Cg = 1/ a/3 in our conformal 
liquid: these two effects to certain degree cancel each 
other. 

Let us further note, that at the initial time the pres- 
sure and flow gradients are especially large at local den- 
sity fluctuations. Therefore the applicability conditions 
of even (viscous) hydrodynamics itself should be inves- 
tigated. Interesting effects, such as e.g. cavitation, are 
known to occur in other hydro applications in similar set- 
tings. Theoretically, the issue is what the sum of all large 
gradients times corresponding dissipative coefficients can 
actually do if resummed: see recent discussion in Ref. [? 
] on "resummed hydrodynamics" and its applicability, 
for AdS/CFT and heavy ion collisions. 

Returning to hydrodynamics, we would like to address 
the issues in the case with maximal symmetry: therefore 
we only discuss the central collisions, which are axially 
symmetric (without perturbations). We also would like 
to be as transparent as possible, thus using the analytic 
tools. Finally, we believe that in this problem, as in many 
others, one should look for the Green function, the solu- 
tion with an elementary delta-function-like source. Once 
it is found, any type of initial conditions can be easily 
included by just a convolution with the Green function. 
From the physics point of view it seems to be more impor- 
tant to calculate the effect of the viscosity on the shape 
of the angular response, rather than to include the non- 
linear interactions between the harmonics, as "event-by- 
event hydrodynamics" does. 

The paper is organized as follows. We start by dis- 
cussing two approaches which can be used in the case 
when the perturbation size is much smaller than the size 
of the system, so that the number of excited harmon- 
ics is large. One is the general "geometric acoustics" 
method, which substitutes the wave equation by mechan- 
ical equations for the propagating "phonons" . The other 
uses the standard eikonal representation of the solution, 
plus separation of variables. Finally, we focus on the so 
called Gubser flow, the conformal analytic solution for 
the fireball expansion [15J with longitudinal and trans- 
verse flows. Significant further development is due to 
Gubser and Yarom [TB] , who derived the linearized equa- 
tions for the propagation of small perturbations around 
it. In our paper we extend their results to sound Green 
function, the coherent sum of all harmonics describing 
the propagating sound from a point- like "hot spot" . Our 
next step is to focus on how the perturbations modify the 
freezeout surfaces and thus observed spectra and corre- 
lators, and finally compare the latter to the data. 



tions of the "average fireball shape" are not included. 
Additional simplifications often used are due to (approx- 
imate) symmetries which the problem possesses, for ex- 
ample rapidity-independence and also consideration of 
only central (axially symmetric) collisions. If those are 
assumed, the number of variables is reduced from 4 to 2, 
and one may start thinking about its analytic treatment. 
Otherwise, the problem only allows for numerical solu- 
tions, which are widely used in practice but will not be 
discussed in this work. 

Our main goal, as we proceed, will be to go to the 
"'first approximation" , deriving small perturbations of 
the zeroth-order solution. Unlike the zeroth order, the 
perturbations are not assumed to have any a priori sym- 
metries. The main object of the hydrodynamical descrip- 
tion, the stress tensor, is conserved: thus the equations 
to be solved are written as its zero covariant divergence 



(I.l) 



where the zero and one in parenthesis are not the indices 
but the order of perturbation. The perturbation term 
will be assumed to be small and treated in the linear 
approximation. 

While it is all very generic, for completeness of the pa- 
per let us remind some details here, starting with the sim- 
plest example of rapidity-independent "Bjorken" flow. 
Even in this case, one needs curved coordinates with 
a non-trivial metric, thus covariant derivatives and the 
nonzero Christoffel symbols will be needed: 
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(1.2) 



Changing Minkowski coordinates t, x, y, z, with z along 
the beam, to the hyperbolic-cylindrical set r, 77, r, 



t = T cosh rj z = T sinh ry 
X — r cos 0, y — r sin 

T = \/ — z^ 7? = - In I 

' 2 V^--^ 

one finds the following metric tensor 
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(1.3) 
(1.4) 

(1.5) 



(1.6) 



and the Christoffel symbols, following from standard ex- 
pression 

F^^. = (1/2)/^ (5,,.,, + 5,,,, - fc) . (1.7) 



A. Relativistic hydrodynamics, the zeroth-order 

By the "zeroth order" hydrodynamical evolution of the 
system we mean the one in which all possible perturba- 



have the following non-vanishing components 

I 



F^ = r 

7]!] ' 
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(1.8) 
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Those are inserted into (1.2) together with the general 



expression for relativistic Navier-Stokes stress tensor 
Tfj^i^ = {e + p)ui_,u^ +P9iiu - 2riaf_,^ - C(u'*';a)A^^, (1.9) 



where, 



U/3-a + Ua-p Qaf} \ 
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(1.10) 
(1.11) 



The first two terms of the stress-energy tensor correspond 
to "ideal hydrodynamics", while the third and fourth 
ones are due to shear and bulk viscosity, respectively. 

The corresponding analytic solution, known as the 
Bjorken flow, [17 corresponds to colliding objects be- 
ing infinite walls of matter, eliminating the transverse 
flow and any dependence on the two transverse coordi- 
nates x,y or r, 0, as well as on rj. Furthermore, we will 
consider the simplest co-moving flow case, with a trivial 
4-velocity = (1,0,0,0). Then the non-viscous stress 
tensor returns to its generic form 

T''^ = dtagieiT),p{T),p{T),p{T)) (1.12) 

in the medium rest frame, depending on the proper time 
T. The resulting 00 and 11 equations, together with the 
thermodynamic identity relating the differentials of these 
quantities 
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e+ p 



(1.13) 



can be put into the final form of one single "entropy 
production equation" 



de 



ds 



dr e + p dr t 



ie+p)T J 



Note that if both viscosities are zero, the solution is just 
ST = const, which implies simply the total entropy con- 
servation. 



II. SOUND PROPAGATION IN THE 
SHORT- WAVELENGTH APPROXIMATION 

A. The geometric acoustics 

If the wavelength of the perturbation is small com- 
pared to the size of the system, one can describe sound 
propagation in the "geometric acoustics" approximation, 
see textbooks such as [THI- The reason we can use such 
an approximation in our problem is the assumed locality 
of the initial "hot spots" (and thus the initial width of 
the propagating circular wave). All we need is that their 
size is much smaller than the fireball dimensions 



I < R 



(2.1) 



The derivation of the approximation is based on the 
analogy between the Hamilton- Jacobi equation for the 
particle propagation and the wave equation for the sound, 
deriving the Hamilton equations of motion for the "sound 
particles" ("phonons" ). The resulting equations of mo- 
tion for them are 



dr 
dt 

dk _ 
~di " 



9a; (fc, r) 
dk ' 

9a; (fc, r) 
df 



(2.2) 



(2.3) 



driven by the (position dependent) dispersion relation 
a;(fc, r). 

Let us start with the simplest non-relativistic case, 
with small velocity of the flow, u <C 1. In this case the 
dispersion relation is obtained from that in the fluid at 
rest by a local Galilean transformation, so that for flow 

u{r) 



U!{k, f) — Cgk + {k u{r)) . 



(2.4) 



As a simple yet relevant example, let us use the (gener- 
alized) Hubble flow in which the velocity profile is linear 



(2.5) 



with some constant (time and coordinate independent) 
Hubble tensor. The eqn (2.31 now reads as "rotation" of 
the phonon momentum 



dr_ 

~dt 



(2.6) 



If the Hubble tensor is symmetric, it can be diagonal- 
ized with 3 real eigenvalues. Hi, H2, H3, so the gen- 
eral solution in its eigenframe is the exponential change 
of the corresponding momentum components ki(t) = 
exp{~Hit)ki{0). Note that if all three eigenvalues are the 
same, the unit vector of the direction would be time- 
independent. Furthermore, if the Hubble tensor contains 
an anti-symmetric part, the direction vector would be 
rotating around the vector CijkHjk- 



Let us now come to the first eqn (2.2) 



dr' 
lit 



= c,nl{t)+W^r^{t). 



(2.7) 



with the first term in the r.h.s. containing a unit vector 
along k. The simplest case is when the Hubble matrix 
is proportional to the unit matrix and the first term is 
time-independent: then the solution is simply a linear 
addition of the sound motion and the Hubble expansion 



r(t) = Cstfij: + r{0)exp{+Ht) . 



(2.8) 



This approximation is enough to explain the deforma- 
tions which the zeroth-order flow induces on the basic 
geometric shapes of the sound fronts - the cylinders, 
spheres or cones - appearing in a non-floating medium. 
(We will use it for this purpose elsewhere [H].) It is 
however not so useful for predicting the corresponding 
amplitudes of the wave, which we will discuss in the next 
subsection. 
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B. Wave equations with separable variables 

Let us explain the idea in the simplest setting, assum- 
ing that there are only time and one relevant space co- 
ordinate, X. Let us also assume that one can eliminate 
the velocity and write the hydrodynamic equations as a 
closed second-order linear equation for the temperature 
perturbation 6{t,x) 



or 

+C3{t,x)— + Ciit,x)S = 



(2.9) 



where C1..C4 are some functions. 

The idea is similar to the semiclassical approximation 
in quantum mechanics, which uses for the wave function 
a form il}{t,x) ~ A{t,x)exp{iF{t,x)/h), with some am- 
plitude and the phase, assuming that the phase is para- 
metrically large F/h^ 1. If so, one can find a solution 
satisfying subsequently parts of the equation of the same 
magnitude. 

Let us show how it works for the generic 2-d equation at 
hand. One also introduces the amplitude and the phase 



S(t,x) ^ A{t,x)exp{i(l){t,x)/e) 



(2.10) 



with the h substituted by a dimcnsionless abstract small 
parameter e. Its substitution into the equation above 
yields three types of terms 



1 



-02 + Ci(0')')] 



i , A(h 
2 — 



A!~ 



A 



C\(\>" -I- C*30'] 



.A A" 



A 
A 



C3^+C4]=0 (2.11) 



For small e one starts from the first square bracket. If 
the first coefficient can be factored into functions of both 
variables, Ci = Cit{t)Cix{x) it can readily be solved 
yielding 



can be only dependent on t, provided C2 depends on time 
only. The last three 0(l/e) terms would be factorizable 
into Cit{t) times a function of a; if C3 = Cit{t) * C^^ix)- 
If so, the solution for both parts of the amplitudes are 



At{t) = exp 

Jo 



Citih) 

4Cit(U) 



-C2t(ti)/2] 



Ax{x) = exp / dxi 



a ^ C[,ix,)) 
y/CiM 4Ci,(xi) 

2Cix{xiy 



(2.13) 



A new separation-of- variable constant a formally appears 
here, but it does not generate anything new in respect to 
what was already included in the phase, so it can safely 
be put to zero. 

Familiar examples of waves are e.g. the spherical and 
conical waves, in which case the variables can be sepa- 
rated . Indeed, when the spatial part of the equation is 
d-dimensional Laplacian, one has 



-, C2 — 0, C3 



1 1 



X d 



(2.14) 



and the corresponding amplitude decays with distance as 



A 



1 



(2.15) 



which is well known for spherical (d=3) and cylindrical ( 
d=2 ) waves. 

As the reader will see later, the sound on top of Gub- 
ser's flow can also be shown to have an amplitude depend- 
ing on new variables p, in a factorizable way, which was 
not the case in the original coordinates, the proper time 
r and r. Therefore, without introduction of these coor- 
dinates, one would not be able to solve the equation for 
the amplitude in such a simple factorized form. 



III. PERTURBATIONS ON TOP OF THE 
GUBSER FLOW 

A. Summary of the Gubser flow 



<l>{t,x) = k[ / ^Cu{ti)dh± 



dxi 



\/Cix{xi) 



(2.12) 



where the separation of variables constant fc, the "wave 
vector", is assumed to be large. When Cix = ^\JC\t = 
Cs = const we have a function of x — c^i, the usual prop- 
agating wave. 

The amplitude A should be found from the second ap- 
proximation, the terms of the order 1/e. One may again 
get an explicit solution assuming the variables can be sep- 
arated. Looking for the amplitude in a factorizable form 
A = At{t)Ax{x) one can see that the first three terms 



The Gubser flow [121 [H] is a solution which keeps the 
boost-invariance and the axial symmetry in the trans- 
verse plane of the Bjorken flow, but replaces the trans- 
lational invariance in the transverse plane by symmetry 
under a special conformal transformation. Therefore, the 
matter is required to be conformal, with the EOS 



e = 3pr^T'^ 



(3.1) 



and the speed of sound = l/\/3- The solution has one 
dimensional parameter q via which the finite size of the 
nuclei is introduced. 
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Working in the (r, 77, r, 0) coordinates with the metric 

ds^ = ~dT^ +T^dr/^ +dr^ +r^d(f>^, (3.2) 

and assuming no dependence on the rapidity 77 and az- 
imuthal angle (jj, the 4-velocity can be parameterized by 
only one function 

Ufj_ = (— cosh k(t, r), 0, sinh k(t, r), 0) (3-3) 

Omitting the details from 15J, the solution for the ve- 
locity and the energy density is 



vj_ ~ tanh k(t, r) 



1 + 9^x2 + g2^2 



(3.4) 



B. Perturbations of the Gubser flow 

Small perturbations to the Gubser flow obey linearized 
equations which have also been derived in [16] . We start 
with the zero viscosity case, so that the background tem- 
perature (no w to b e called T;,) will be given by just the 
first term in (3.11). The perturbations over the previous 
solution are 



efined by 



T 



7^(1 + ^) 



(3.12) 
(3.13) 



with 



(2g)8/3 



r4/3 (1 + 2(72(r2 + r2) + g4(^2 _ ^2)2)4/3 



(3.5) 



where Eq is some normalization parameter. 

In |16) Gubser and Yarom re-derived the same solution 
by going into the co-moving frame. In order to do so they 
rescaled the metric 



ds^ = T^ds^ 



(3.6) 



and performed a coordinate transformation from the r, r 
to a new set p, given by: 



1 — g^T^ + q^r'^ 
2qr 



sinh p = — 
tan 6 = 



1 + (727-2 _ q2j,2 
In the new coordinates the rescaled metric reads 



(3.7) 
(3.8) 



ds-" 



'dp^ + cosh^ p (d0^ + sin^ ed(f>^) + dTf (3.9) 



and we will use p as the "new time" coordinate and 9 as 
a new "radial" coordinate. In the new coordinates the 
fluid is at rest, so the velocity field has only nonzero Up. 

The relation between the velocity in Minkowski space 
in the (r, r, (/), "q) coordinates and the one in the rescaled 
metric in (p, 9, 4>, rf) coordinates corresponds to: 



(3.10) 



while the energy density transforms as: e = T~^e. 

The temperature (in the rescaled frame, T — rfl^^T, 
with = e/T'^ = 11 as in [TS]) is now dependent only 
on the new time p, and in the case with nonzero viscosity 
the solution is 

. to , Hosmh3p „ /3 7 5 „ ^ 

T = H -7- nFi — Sinn p 

{coshp)2/3 9(coshp)2/3 6 2' 



(3.11) 

where Hq is a dimensionless constant made out of the 
shear viscosity and the temperature, i] — HqT^ and 2F1 
is the hypergeometric function. In the inviscid case the 
solution is just the first term of expression (3.111, and 



of course it also conserves the entropy in this case. The 
picture of the explosion is obtained by transforming this 
expression back to the r, r coordinates and performing 
the appropriate rescaling. 



(-1,0,0,0) (3.14) 
(O,7/e(p,0,0),u^(p,0,0),O) (3.15) 
5{p,9,cj)) (3.16) 



The careful reader will notice here, that although general 
perturbations should not have any symmetries of the ze- 
roth solution, we have not listed rapidity among the vari- 
ables. Indeed, we only consider the perturbations which 
are rapidity-independent. The reason for that is that the 
initial state perturbations are initiated in the transverse 
plane but rapidity-independent, so that the waves they 
induce also propagate in the transverse plane only. 

Plugging expressions ( |3.12 |,(3.13) into the hydrody- 
namic equations and only keeping linear terms in the 
perturbation, one can get a system of coupled 1-st order 
differential equations. Furthermore, if one ignores the 
viscosity terms, one may exclude velocity and get the 
following (second order) closed equation for the temper- 
ature perturbation: 



dp^ 
4 



1 



3 cosh^ , 



(Pi 

96*2 ^ tan( 



1 85 



1 8^5 



89 sin^ 86^ 



1 95 n 
tanhp— = 

3 8p 



(3.17) 



As we will show, it has a number of remarkable proper- 
ties. 



C. The short-wavelength approximation for the 
sound waves on top of the Gubser flow 

Before we proceed to the exact solution of this equa- 
tion, let us follow the procedure described in section IIB 



and study the solution to equation (3.171 in the short 
wavelength approximation. We start by looking for a 
factorized solution of the form: 

S = e'^^^^P^-f'^'^^-^*'^^^^Fp{p)Fg{9)F^{(j>) (3.18) 

where fi >> 1, such that the derivatives taken over the 
exponential are dominant. In this way, we study the 
equation separating it in different equations depending 
on which power of the derivatives over the exponent they 



have. The first step is to look only at the second deriva- 
tives because, since they produce terms of second order 
in the exponent, they are the leading ones. In this way 
we find: 



fpip) 
feiO) 



± —pz k arctan e'' + A 



= ± / ddx k'^ - 



B 



sm 



= ±7710 + C 



(3.19) 

(3.20) 
(3.21) 



The integral in (3.20) can be solved, but it gives a cum- 
bersome result. So in what follows (of this section) we 
will assume no <f> dependence just to get an idea of the 
result. When we do this, the functions in the exponent 
reduce to: 



foip) — i^^fcarctane'' 
fg{e) = ±ke + B 



A 



(3.22) 

(3.23) 
(3.24) 



The function fp{p) is almost linear in p in the region 
that we are interested in studying (—2 ^ p ^ 1), so we 
find the phase of the solution to be ^ kp which means 
that we indeed expect to find solutions in the form of the 
sound wave propagation (in this region). 

Now that we have found the functions in the expo- 
nent we look for the wave amplitude by cancelling among 
themselves the terms with the first power of the large ex- 
ponent: by doing this we find the amplitudes to be 



Fe{e) 



1 



(cosh p)i/6 

1 



(3.25) 
(3.26) 



D. 



The exact separation of variables for the 
perturbation 



We have seen that in the short wavelength approxi- 
mation we found a separable wave-like solution to equa- 



tion (3.17), and now we would like to see if the ex- 



act solution can be found by using variable separation 
S{p,9,(j>) = R{p)e{0)<i>{e). It is indeed so. In the non- 
viscous case, that we are now discussing, each of the three 
equations 



R"{p) + ^ tanhpi?'(p) + ^ Rip) = (3.27) 
o 3 cosh p 



^ ' tan6' ^ ' \ sin^e 



] 6(61) = (3.28) 



$'(0) -l-77i2$((/)) = (3.29) 



is analytically solvable, with the result 



R{p) = 



2/3 



2/3 



-(tanhp) + C2Q"'i 



+ iVl2A+ 



-(tanh p) 



(cosh pY^^ 
0(0) = CsPr {cos 9) + C4Qr {cos 9) 



(3.30) 



where A = Z(/ -|- 1) and P and Q are associated Legendre 
polynomials. The part of the solution depending on 9 
and (p can be combined in order to form spherical har- 
monics Yi„i{9,(l)), such that 5{p,0,(l)) cx Ri{p)Yim{0 , (j)) . 
This property should have been anticipated, as one of 
the main ideas of Gubser has been to introduce a coordi- 
nate which together with (p make a map on a 2-d sphere. 

The implications of that for the physics we are going to 
discuss are as follows. While we will project the spectra 
and correlation function to the azimuthal angle 4> and its 
Fourier components, we will be focussing on the quantum 
number m conjugated to it. In particular, the community 
is very much focused on the "triangular flow" with 771 = 
3. In principle, however, this is produced by many l- 
harmonics, providing the obvious condition / > 771 holds. 
Harmonics with different / have obviously different radial 
dependence. (We mention this point, because there was 
some controversy about the powers of r, especially for 
various definitions of the "dipole flows" with 7n = 1.) 

Let us explore the asymptotic behavior of the Legendre 
functions when Z >> 1 that is given by [5D]: 



P/"(cos^ 
OJ"(cos^ 
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r r(/ 
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niTT > 
2 ' 



rriTT > 
2 ' 



r(/ + 3/2) 



\/2 sin 9 



(3.31) 



These expressions show that for large I the solution 
presents oscillatory behavior in 9 with an amplitude given 
by ^^r\—g - It gratifying to see, that this is the same 
that we obtained in the short-wavelength approximation 



for Fg{9) (eq 3.261 in the previous section. 

Now let us look into the p-dependent part of the solu- 
tion in the large I limit we have that the Legendre poly- 
nomials as a function of tanhp correspond to: 



(tanh p) 
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2 r(I -I- m + 1) 
r(i + 3/2) 



\/ cosh p 



1 \ 7T TTITT 

I -\ arccos (tanh p) — — + 
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\/ cosh p 



1 \ TT niTT 
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(3.32) 



Again we see an oscillatory behavior and a wave am- 
plitude. In this case the amplitude is given by V cosh p 
and if we divide this by (coshp)^/'^ as we have in the ex- 
act solution (3.301 we get an amplitude for the wave of 
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(cosh pYI<' 



which is the same as we got in the preceding with 



section (3.251 using the short wavelength approximation. 



So we have checked that for large / (5(p, 6*, 0), and there- 
fore the temperature perturbation in the rescaled frame, 
Ti{p,9,(p) = Tb{p)5{p,9,(j)), does behave like a sound 
wave. 



E. Propagation of the local initial-state 
perturbation 

Let us study the propagation of the hydrodynamical 
response induced by an initial perturbation on top of 
the background at some initial "time" po, given by a 
Gaussian-shaped initial "hot spot" : 



fi{pQ,e,(j)) (X e 



(3.33) 



We further assume that at the initial time there is no 
flow (momentum), only extra energy, so another initial 
condition is: 



ue{pa) = 
u^{pa) = 



(3.34) 



which define the initial derivative of the temperature per- 
turbation, since [16] 



vi{p) = 



(3.35) 



(3.36) 



Uli = Vi{p)drYim{0,4)) 

3 cosh^ p dSi 
+ iLp 
where i = 6,(j). Thus we require 

^1 - 

The general solution for linear perturbations is 

711 — I 

flip, 9,^) ^ ClmRl{p)Ylm{9,(l>) (3.37) 

I m— — l 
ni—l 

uix{p,0,(p) = X! X! cimVi{p)diYi^{9,(l}) (3.38) 



I m— — l 



with 



Riip) 



A,P 



.(tanh p) + BiQ''^f 



- 1 + i yi21(i+l)T 



-(tanh p) 



(cosh p)'*/3 



(3.39) 



where c/™, Ai and Bi are constants that can be deter- 
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(tanh p) 



(cosh p)2/3 
^^■^ (tanh p) 



Q 



-; + g\/l2U' + l) + l 

(cosh p)2/3 



(3.42) 
(3.43) 



which, together with the Yim{9,(j)), finally provides the 
complete solution for S in equation (3.171 . The denomi- 



nator of the right term of Ri (p) is the so called Wronskian 
of the functions pi{p) and qi{p) evaluated at the initial 
"time" pq. Since the Legendre polynomials P and Q are 
linearly independent, the Wronskian is always non-zero, 
so we are guaranteed that the function Ri{p) is always 
finite. 




FIG. 1: (Color online) Top: Ri{p) defined in for har- 

monic number I from 1 to 10 (from less to more oscillating 
ones). All curves are arbitrarily normalized to 1 at p = —2.07. 
Bottom: the corresponding harmonics of the velocity vi{p), 
defined in (3.351, also for I from 1 to 10. Note its oscillatory 



mined using the initial conditions (lOsl) and (fOel). With behavior for larger / and later times 



Ai and Bi determined, the /9-dependent part of the tem- 
perature is 

/ cosh Pq ^ 



Riip) 
Slip) 



\ coshp 



Slip) 



tUp^(p)-^lpogKp) 

tloPliPo) 



dp 'Po 



mipo) 



(3.40) 



(3.41) 



The first ten harmonics for Ri and vi are plotted in 
Figjl] showing how the amplitude varies as a function of 
"time" p. One can see how the initial deformation (the 
upper plot, set to one for each 1 for comparison) is trans- 
ferred into the flow velocity (the lower plot). One can 
also see that while for the lower harmonics it happens in 
a more or less linear way, higher harmonics show oscillat- 
ing behaviour , as expected for sound waves. Indeed, one 
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should see a transition from potential to kinetic energy 
happening with higher and higher rate, as the harmonic 
number grows. 

The cim coefficients are calculated using the orthogo- 
nality of the Legendre polynomials, and are given by: 



Clr. 



2tt 



(3.44) 



Once we find all the constants, we can study the evolution 
of the perturbation given by expression (3.371. In figure 
[2] we show three frames from a movie-like evolution in 
T of a perturbation Ti(t, r, (^) produced by a local "hot 
spot" which was calculated for a Gaussian centered in 
6 = 1.5 with a small size s — 0.1, which corresponds to 
a perturbation localized at r = 4.1 fm and with a width 
of 0.4 fm. Notice that while the perturbation is in the 
rescaled frame, we are using the regular coordinates r, r. 

We have used 30 harmonics for this movie, and it is 
nice to see that they all add up coherently into a con- 
sistent picture of a sound wave propagation. While it 
does correspond to a qualitative picture of a circle from 
a stone thrown into the pond, with which we had started 
this work, it is in fact an exact solution, riding on the 
zeroth order explosion picture which is by itself rather 
complex. In order to find the analytical expressions for 
the perturbation on top of the fireball it was necessary to 
invent the p and 9 coordinates, so that all of the expres- 
sions can be factorized in terms of these coordinates. So 
a lot of correct thinking was needed, to make this movie 
possible. 



F. The viscous effects 

In the second paper of this series [T| we introduced the 
viscosity-based scale, which all structures produced by 
point-like perturbations would obtain at freezeout. With- 
out going into details, let us just remind the reader that 
while the width of the circle grows with time as t^^^, 
its radius grows as r, and therefore the relative contrast 
(the former divided by the latter) is improving as r~^/^. 
As far as the amplitude of the wave is concerned, in a 
short-plain-wavelength approximation the stress tensor 
harmonics with momentum k are attenuated by a factor 



6T^^{t, k) = exp 



2 T] kH 
'3 



5T^,,(Q,k) (3.45) 



known from textbooks on sound, sometimes called "the 
viscous filter" . Note that its exponent contains the mo- 
mentum squared, due to the extra derivative in the vis- 
cous tensor, and therefore the effect of viscosity for the 
higher harmonics is strongly enhanced. Obviously, the 
same qualitative behavior is expected for our Z,m har- 
monics. 

The basic equations for the p-dependent part of the 
perturbation, now with viscosity terms, can be written 



r,o.5 




r,o.5i 
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FIG. 2: (Color online) Evolution of the perturbation in the 
rescaled frame but in the regular coordinates Ti (r, r, <j)) from 
(3.371 using the change of coordinates 3.7 and 3.8 From top 
to bottom:r = 1 fm/c ,t — 4 fm/c ,t — 6 fm/c 



as a system of coupled first-order equations [T6j. We are 
assuming rapidity independence, thus the system of equa- 
tions (107), (108) and (109), from the referred paper, be- 
comes two coupled equations, for (the p-dependent part 
of) the temperature and velocity perturbations 



dw 
dp 



= —Tw , w = 



(3.46) 



where the index v stands for viscous and the matrix com- 
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ponents are, 
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(3.47) 



+ 6Hi tanh^ / 



(^Ho tanh p 



2T,, 



This system can also be written as a closed second order 
differential equation for Sy (p) : 

5v dSv 
dp'^ dp 

Til dri2 



rii 



1 rfFia 

ri2 dp 



(3.48) 



fdTi 

V dp ri2 dp 







Unfortunately, unlike the zero viscosity case considered 
above, the equations one gets after separation of vari- 
ables cannot all be solved analytically and thus have to 
be solved numerically which has been done using Math- 
ematica's ODE solver. The part of the solution which 
depends on 9 and (j) is not affected by viscosity, so it con- 
tinues to be given by the spherical harmonics Yim{9,(j)). 

Our results for the nonzero viscosity will use either 
Hq = 0.33 (7?/s = 0.134), like in [T6j, or the value = 
0.19 {rj/s = l/(47r) = 0.08), the conjectured lowest value 
possible predicted by AdS/CFT in the strong coupling 
limit. 

In Figj3]we plot the "time" dependence Syi{p) for sev- 
eral harmonics and compare them to the inviscid case 
Si{p) for some I's. As expected, the viscosity reduces 
higher harmonics more, but as far as time dependence 
is compared to inviscid case, we see that viscosity liter- 
ally kills the contribution at certain time, which becomes 
shorter and shorter for larger I. As the time is limited 
by the freezeout time, we observe that the contributions 
of all sufficiently large I > Imax ~ 10 become completely 
negligible. 

The p-dependent part of the velocity can be calculated 
once Sy is known: 



Vvl{p) 



51 



vl 



(3.49) 



For the first 10 1, the curves Vvi{p) are plotted in Figj4] 
Comparing this to vi (p) at zero viscosity (bottom plot of 
Figjl]) we see that the amplitude for the velocity is also 
damped in the viscous case for large 1 and increasing p. 



IV. APPLICATIONS 



A. Matching the Gubser flow with the heavy ion 
colhsions 



With the exact solution to the perturbation equation 
riding on top of the Gubser flow at hand, one may go 
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FIG. 3: (Color online) Comparison between the magnitude of 
harmonics Si{p) (3.41 1 for the ideal case (soUd blue lines) and 



Svi{p) (3.481 for the viscous case with rj/s = 0.134 (dashed 
red lines), for I = 1,3, 5, 10, from top to bottom. 



back to the r and r coordinates and try to calculate what 
should happen in real heavy ion collisions. 

But before we do so, let us remind the reader once 
again that the Gubser flow is by itself an idealization of 
reality. The real hadronic matter can only be approx- 
imated by the conformal EOS e — 3p during its QGP 
phase, which lasts about 1/3 of the total time at RHIC 
and perhaps around 1/2 time at LHC. The rest is the 
near-Tc domain and the hadronic phase, in which the 
speed of sound changes from 1/a/3 = .577 to about .35 
and .45, respectively. Although this change is not very 
large, we do notice that the radial flow obtained with 
the Gubser flow is too large. Respectively, the freezeout 
time Tpo is indeed somewhat smaller than that observed 
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FIG. 4: (Color online) The velocity harmonics i (p) (defined 
in (3.49 1) for the first ten values of /, in the viscous case with 
77/s = 0.134. To be compared with figure [T] bottom. 



in numerical hydrodynamics with correct EOS. Perhaps 
some of our results for the perturbation should also need 
some adjustment, due to these facts. 

The second similar comment is that Gubser's solution 
has a particular shape, which has no reason to coincide 
with the shape of the real Au nuclei. The finite size of 
the fireball and the shape of its temperature profile is 
determined by the parameter q which we take equal to 
(4.3/m)~^ following [TS]. The second parameter that we 
need to fix is the constant Tq . Again from [T3] we get the 
formula for this parameter 



with 



Tn 



e 



1 



3 dS 



11, 



dS 
dr] 



7.5- 



dNr. 



drj 



(4.1) 



(4.2) 



For central (0 - 5%) coUisions at LHC dN^h/dr] ^ 1600 
PT] which gives a value of Tq « 7.3. Using these val- 
ues one gets a freeze-out time t/o ^ 6, which is rather a 
short time that doesn't allow for the sufficient evolution 
of the sound circles. Since we are interested in study- 
ing the propagation of sound perturbations and the size 
of the sound horizon depends on the freeze-out time, we 
will use Tq « 10.1, which corresponds to having about 
2.Q{dNch/drf)ijHC- It is important to note that the back- 
ground temperature in the ideal case corresponds to 



T 



1 



To 



(coshp)2/3 ' 



(4.3) 



so we are using an initial temperature of about 630 MeV. 
The parameters we used are such that the size of the fire- 
ball at freeze-out, the radius of the sound circle and over- 
all transverse expansion velocities v_\_{r,t « tf) mimic 
reality of RHIC/LHC collisions. The price for that is 
somewhat too large initial temperature and overall en- 
tropy. 



The hydrodynamical equations should be used only af- 
ter some approximate equilibration of hadronic matter is 
achieved. While the mechanism of it, as well as precise 
timing remains unknown, we do know its order of mag- 
nitude to be a fraction of fm/c. For our calculations 
we assume that thermalization occurs at the initial time 
T = 1 /m/c, and it is at this time that we define our ini- 
tial "hot spot" and start evolving it using hydrodynam- 
ics. One can do so until the final freeze-out is reached, at 
which point the interaction between secondaries becomes 
ineffective and sound propagation stops. Below we will 
discuss how the hydrodynamical perturbations should be 
translated into the experimental observables. 

Let us point out that we study the effect of a single 
hot-spot on the fireball which we characterize as a Gaus- 
sian temperature perturbation on top of the background 
temperature. In real collisions, there are many such per- 
turbations, but since we solve the problem in the lin- 
ear approximation, their evolution is mutually indepen- 
dent. Furthermore, in the experimental statistical study 
of small two- or three-particle correlations, the contri- 
bution of the uncorrelated fluctuations is cancelled out 
automatically. 

In Figj2]we see that at the time t — 1 fm/c a Gaussian 
"hot spot" centered at r = 4.13/m, (j) = tt corresponds to 
having it at "time" p — —2.07 centered aX 6 = 1.5, = tt. 
Of course, since p — p{t, r) at any given time t, p depends 
on r, so the initial condition p = —2.07 is for the center 
of the Gaussian. 



B. Modification of the freezeout surface and of the 
particle spectra 



The standard expression for a spectrum. 
Cooper- Frye formula [22], is given by 



E 



,dN 
d^p 



= - / 



T 



known as 



(4.4) 



The overall minus is there because we work using the 
mostly plus signature. The function f corresponds to 
the thermal distribution inside the fluid cells, boosted by 
their hydrodynamical motion at the time of the freezeout 



fip) 



1 



exp(—pf^u^/T) ± 1 



(4.5) 



for Bose/Fermi particles. (In reality, we will be only in- 
terested in the tail, so the Boltzmann approximation will 
always be enough.) The minus sign in the exponent is 
because we are working in the mostly plus signature. 

The temperature and velocity in this formula are sup- 
posed to have a space-time dependence derived from hy- 
drodynamics. The freeze-out surface E'' that appears 
in the Cooper-Frye formula corresponds to a certain ki- 
netic condition, of the form that the ratio of a particular 
reaction rate to the matter expansion rate reaches a par- 
ticular value. Since there are many reactions involved in 
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the process, strictly speaking there are multiple freezeout 
surfaces. One usually separates "chemical" and "kinetic" 
freezeouts, in which inelastic and elastic scattering rates 
are involved. Since different secondaries (pious, K, uu- 
cleons, ... J/tp) in fact have quite different elastic cross 
sections, the "kinetic" surfaces should in fact be different 
for each species. 

We are not going to discuss all those complications in 
this work, and think of only one type of secondaries, the 
pious. Furthermore, we will use a drastic simplification 
often used, assuming that the freezeout surface is the 
isotherm T(t,x) — Tpo- If so, the surface can be de- 
termined from hydrodynamical output, for example its 
time-like part can be written as 



{Tfoix,y),x,y,r)) 



(4.6) 



where t/o is the time at which the fireball reaches the 
freeze-out temperature. The Cooper-Fry formula con- 
tains the vector normal to the surface which is then 



dT," dT,^ d^P 



dxdydrj (4.7) 



dx dy drj 

Tfodxdydi] (4.8) 



dx dy 

Here g is the determinant of the metric and e^^^Ap is the 
Levi-Civita symbol. 

The perturbations affect the spectra in two ways. 
First, the flow velocity in the exponent is corrected by the 
extra terms of the first order due to sound. The second 
effect, related with the first order temperature perturba- 
tions (l + S), are more subtle. Hotter matter (positive S) 
in the event with a "hot spot" and perturbation from it 
imply a production of extra entropy density (increases by 
(1 -I- S)^) as compared to the zeroth order fireball. This 
means there would be extra secondaries produced, as this 
entropy is "hadronized" . By assumption, it happens lo- 
cally, delaying a bit the freezeout according to condition 



To{t,x) [l + 6{x,t)]^T, 



FO 



(4.9) 



Thus delay is absolutely necessary, it provides extra vol- 
ume for the extra matter produced, as compared to the 
zeroth order explosion, since by assumption the freeze- 
out temperature and thus the matter density at the FO 
surface are held constant. The deformation of the FO 
surface not only increases the volume, giving place for 
the extra particles just discussed, but it also prolongs 
hydro evolution, providing a bit larger flow. 

Let us now discuss another issue: at what part of the 
particle spectra we should focus, in order to see best the 
effect of the perturbation. The Cooper-Fry formula has 
Pt of the particle in the exponent, so it is tempting to 
take it as large as possible. And indeed, all hydro effects 
(such as e.g. the elliptic or radial flow) are enhanced by 
the increase in the particle momentum pt . There are two 
practical limits to an increase in pt, however: 
(i) One can be understood inside the hydrodynamics it- 
self. The viscous term has an extra gradient, relative to 



the ideal part of the stress tensor. This means that the 
relative role of viscous corrections will grow with pt, till 
at some point it will no longer be small as compared to 
ideal term. Obviously at such pt hydrodynamics should 
be substituted by some other tool, e.g. some kinetic the- 
ory description. 

(ii) In real collisions some secondaries originate from hard 
scattering and subsequent jets. In spite of significant jet 
quenching, at large enough pt the hard component of the 
spectra supersedes the hydrodynamical spectra. Obvi- 
ously, beyond this point one looses ability to follow the 
hydrodynamical component. 

The transition between the hydrodynamic part of the 
spectrum and the hard QCD tail has been determined 
to be between 4-5 GeV [221 HI] so, a bit conservatively, 
we will consider pt = 1 GeV, as a region well inside the 
hydrodynamical domain. Even at this pt, its ratio to the 
kinetic FO temperature is a large number Pt/Tf = 0(10), 
which can be treated as a large parameter of the problem, 
residing in the exponent. 

Let us work out the first-order corrections appearing 
from the perturbation. There are two effects, one from 
the extra matter T = Tf -\- 5T and one from extra motion 
of the matter in the sound wave. The latter contribution 
comes simply from adding the perturbation to the veloc- 
ity, 



5u„ 



(4.10) 



(5m^ is the perturbation, written in (3.38) as ui 



times r. 

The effect due to the extra matter is included when 
calculating the freeze-out surface: 



Tfo = Tf,(T,r) +(5T(T,r, ( 



(4.11) 



where 5T 



Ti/t, with Ti from ( |3.37| .The equation 
(4.11 ) is solved for T(r, 0), and the result for the inviscid 



case IS presented in FigjSj Since the contribution from the 
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FIG. 5: (Color online) Freeze-out surface t{x, y) for the in- 
viscid case. 

perturbation is small, we write r(r, (j)) = Ti,{r) + Srir, (p) 
and consider terms up to first order in ST{r, (p). By this 
we mean that the exponent will be approximated by 

P^u^jn + St) p'^ubf.in) 1 dipf^ub^jn + St)) 

Tf " Tf +Tf d{5T) '^^=°'" 



p''SUf,(Tt) 

Tf 



(4.12) 
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Fig|6] shows St for both, the inviscid and for the viscous 
case. In the former case the contribution is much larger 
than in the latter, where the viscosity has damped and 
widened the peaks. 






FIG. 6: (Color online) Excess of freeze-out surface ST(r,(f)) 
due to the initial perturbation. Top: ideal case, bottom: vis- 
cous case with "q/s = 0.134. Only the half of the surface that 
is affected by the presence of the perturbation was plotted. 



Figure [7] compares the particle distribution for three 
cases, (i) the inviscid case, (ii) the minimal viscosity case 
•q/s = l/(47r) and (iii) the case where r]/s = 0.134. In 
the ideal hydro case the two peaks of the angular dis- 
tributions, due to the overlap of the perturbation with 
the fireball boundary, are more pronounced than in the 
cases with nonzero viscosity. Also, in this case (i) one 
can clearly see high frequency oscillations on the curve. 
Those are an artifact of the arbitrary limit of the number 
of harmonics used to I < Imax = 30. The oscillations dis- 
appear when we take viscosity into account, because, as 
we mentioned earlier, viscosity kills all higher harmonics 
anyway, with I > Imax ~ 10. In the presence of viscosity, 
the peaks in the particle distribution are weakened, and 
their angular separation is a bit more spread than in the 
inviscid case. 

C. Two-particle correlations 

Looking at experimental data on normalized two parti- 
cle correlations, such as the one shown in the last plot of 
Fig|8j one sees that the peaks are of the order of about a 
percent. This means that the perturbations to the back- 
ground are small, and such small changes cannot be ob- 
served on an event-by-event basis, but only in a large 
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FIG. 7: (Color online) Single-pion distribution in arbitrary 
units as a function of the azimuthal angle (rad), for trans- 
verse momentum pr ~ 1 GeV and rapidity y = 0. From top 
to bottom, the curves are for different viscosity-to-entropy 
ratios, ry/s = 0, l/{4n), 0.134 respectively. 



sample of events. This is why the observables are the 
two(or more)-particle correlation functions, in which the 
non-trivial correlations are separated from the uncorre- 
lated background. Note, that not only fluctuations in 
different events are uncorrelated, but also statistically in- 
dependent fluctuations at different locations in the trans- 
verse plane in the same event. 

In the two-particle correlation functions one measures 
mean squares of the perturbations. Therefore the small- 
ness of the perturbation appears quadratically, and thus 
one has to be able to get to the level between 10^"^ and 
10~^ or so in the correlation magnitude. Nevertheless, 
the large set of the recorded events (~ 10^) by RHIC 
or LHC detector, with ~ 10'^ particles or ^ 10^ parti- 
cle pairs per event provides a sufficient statistical data 



14 



sample. 

Let us now proceed with our theoretical calculation of 
the two-body correlation function based on the single- 
particle distribution resulting from the Green function 
(point-like perturbation) . These correlation functions are 
presented in two forms, which in fact contain equivalent 
information: as a function of the relative azimuthal angle 
or as a "power spectrum" of the flow harmonics. Let us 
start by looking at the former. 

In order to calculate the two-particle distribution one 
should simply take a product of two single-particle dis- 
tributions, and perform the averaging over the random 
axial position of the initial perturbation 

/" ill 

The averaging reduces the function of two angles into 
a function of only one, the azimuthal difference A0 = 
01 — 4'2- (This is only so for central collisions, which are 
axially symmetric: otherwise the situation is more com- 
plicated as the direction of the impact parameter breaks 
the axial symmetry. This is one of the reasons we focus 
on central collisions in this work.) 

Our results for the two-particle distributions for three 
viscosity values are shown in the top three plots of Fig|8] 
Note first their distinctive shape, with a larger peak cen- 
tered at A(j) = (when both particles belong to the 
same maximum of a single-particle distribution) and two 
smaller peaks at Acf) ^ ±2, when two particles belong 
to two different peaks, connected by some flat region be- 
tween them. This shape of the sound Green function is 
in fact very similar to what is observed experimentally, 
for example in the bottom plot in |8] which corresponds 
to data from ATLAS [15]. 

Now comparing the three pictures in more detail, one 
observes that the upper plot (for zero viscosity) has more 
structure. The upper plot has four distinct "dips" in 
which that two-particle distribution is less than average. 
Their origin is explained by matter sucked out by the 
passing sound front behind it. 

The origin of the additional peaks next to the zero- 
angle one is the correlation between one of the peaks 
in the single-particle distribution with matter inside the 
circle. These extra peaks are attenuated when viscosity is 
used and for rj/s — 0.134 they have already disappeared. 
This happens because the viscosity induces cancellations, 
between the negative "suction regions" and positive extra 
matter inside the circle. 

There are now many experimental results for the two 
particle correlations in central collisions such as STAR 
collaboration data [25] for a centrality of — 12%, data 
from ATLAS and ALICE in the very central region 0—1% 
[IS] , [57] . Now, comparing our calculated two-particle dis- 
tributions Fig [8] to these data one should be impressed by 
a striking similarity between their shapes, especially for 
the "realistic viscosity" (the third in Figjs] ) . The width 
of the main peak is correctly reproduced, provided the 
viscosity is correct. Also the "double-hump" structure 
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FIG. 8: (Color online) The two-pion distribution in arbitrary 
units as a function of azimuthal angle difference A(j> (rad). 
From top to bottom, the three upper plots correspond to cal- 
culated shapes, for values of the viscosity-to-entropy ratios, 
rj/s = 0, l/47r, 0.134, respectively. The bottom plot, shown 
for comparison, is measured at LHC by ATLAS collabora- 
tion [25], for the most central (0-1%) collisions. Very similar 
shape of the correlation function is in fact observed by all five 
collaborations at RHIC and LHC. 

on the away side, with the correct shape of the plateau 
in between is found. (The peaks are a bit shifted, it is 
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because the sound velocity as well as the shape/size of 
the freezeout surface is not quite realistic in our analytic 
approach.) 

Let us emphasize that this non-trivial shape conies 
from the hydrodynamical calculation itself, with the ini- 
tial condition simply being a local(delta function like) 
"hot spot" . This agreement of the shape allows us to con- 
clude, that the experiments in question do see the sound 
waves propagation, by a distance comparable to the fire- 
ball radius. The angular positions of the secondary peaks 
depend entirely on the ratio of the "sound horizon" to the 
size of the fireball (the speed of sound and the freezeout 
time). 

All our pictures are assumed to be rapidity indepen- 
dent, thus the zero-angle peak is nothing else but the so 
called "soft ridge" discussed in literature as a separate 
phenomenon. We are pleased to see that its height, with 
respect to the two other peaks, is about the same as in 
the data, especially for the third case in Figjs]). The an- 
gular width of this main peak is, in this case, also quite 
close to the data. 



D. The power spectrum and the initial width of 
the perturbation 

We have also calculated the so called "power spectra" 
for the two-particle correlation functions. Those either 
can be calculated from the Fourier transform of the cor- 
relator as a function of A0, C„, or from the modulus 
squared of the flow harmonics in the single particle spec- 
trum, since C„ = v^. In this last form the expansion of 
the two-particle correlation function is 

^ - l + 2^|v„,|2cos(mA(/.). (4.14) 

and thus it carries the same information as the power 
spectrum of harmonics, in which |v„ip are plotted ver- 
sus m. (Notice that these v„i are the coefficients of the 
Fourier expansion of the particle distribution and are 
not to be confused with the velocity coefficients w/(p) of 
the perturbation). The main advantage of studying the 
power spectrum is that the phenomena associated with 
higher harmonics becomes more visible, which is difficult 
to see in the correlation function itself. 

The result is shown in Figj9] and it presents maxima 
and minima. This structure of the power spectrum, with 
several "acoustic peaks", is known also for other oscil- 
lations, most notably for those seen in the power spec- 
trum of the angular harmonics of the Cosmic Microwave 
Background (CMB) distribution over the sky such as the 
famous Fig. 9 of Both in the Big and Little Bangs, 
the time allocated to the hydrodynamical stage of the 
evolution is limited by the so called "freezeout time" t/, 
after which the collision rates in matter can no longer 
keep up with the system's expansion. At this time the 
propagation of the elastic waves stops and each harmonic 
has at this moment a different phase of its oscillation. 



While the CMB measurements read the temperature 
perturbation S{fo) directly from the sky, and thus the 
nodes of 5i{fo) correspond to the minima, in the Lit- 
tle Bang one has to calculate the specific combination of 
the temperature and flow perturbations. This includes 
the calculation of how the freezeout surface is modi- 
fied, which was done in preceding sections. It is the 
nodes/maxima of this "observable" combination which 
make the acoustic minima/maxima. Note that the sim- 
ple physics behind this argument makes it very robust. 
The minima/maxima are easily predictable and rather 
insensitive to many details such as dissipation. In fact 
the only assumption needed for this idea to be used in 
practice is that the initial state perturbations Si{in) do 
not have an oscillatory dependence on I of their own. 

Before we discuss the results, we need to mention an- 
other important parameter of the problem, namely the 
size of the initial perturbation. In all the discussion 
above this was taken as small and thus unimportant: 
one could think of the perturbation as being practically 
point-like, and thus the results being basically the Green 
function of the equations we are solving. However, as we 
will see shortly, when one discusses the magnitude of the 
higher harmonics, this size does matter. 

FigjO] shows how this works in practice, the three plots 
correspond to three different widths of the initial pertur- 
bation: 0.4,0.7 and 1 fm, and as one can see a change 
in this size does change significantly the tail of higher 
harmonics, the larger the width the smaller the height 
of the larger harmonics in the power spectrum. Never- 
theless, this does not affect the location of the acoustics 
dip and the secondary maximum, which remain around 
m — 7 and 9, respectively. 

Different curves on the plot correspond to different vis- 
cosities (see the caption), and as one can see, they do af- 
fect higher harmonics drastically. This is to be expected, 
as higher harmonics of the flow have higher gradients of 
the flow. One can also see from these figures that the fit 
to the viscosity value must be done together with the fit 
to the initial size, as they are very much correlated with 
each other. 

We will not attempt an actual fit here, adding just 
some comments about the issues encountered. The 
physics of the initial perturbation size should be, first 
of all, related to the size of the "gluonic spot" in a 
nucleon, propagated via pQCD evolution to appropri- 
ate X and scale Q under consideration. At RHIC, with 
X - 10^2 Q ^ 1 - 2GeV we know from DESY exper- 
iments (e.g. diffractive J/ip production) it to be rather 
small, of about .S fm. But then there is some non- 
equilibrium stage, before hydro equations become valid, 
during which this spot should grow. To define the partic- 
ular value one needs to know the non-equilibrium physics 
at this stage. Even to define the start of hydro, one needs 
to know which version of hydro is used, ideal, viscous or 
"resummed" : for recent discussion of these issues refer 
to [? ] and references therein. One more comment on 
the plots in Figj9] is perhaps in order: as the reader can 
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FIG. 9: (Color online) Spectral plots for three for three 
widths of the initial perturbation, 0.4,0.7 and 1 fm, from 
top to bottom. The (magenta) small-dashed, the (red) dash- 
dotted, the (green) solid and (black) dashed curves are for 
r)/s = 0,0.08,0.134,0.16, respectively. The data points are 
preliminary data from ATLAS reported at QM2001 [25]. Sim- 
ilar data (not shown here) have been reported by the PHENIX 
|28| and STAR [29] collaborations. All the curves are arbi- 
trarily normalized to fit the third harmonic. 



see, the curves look shifted toward the larger m from the 
data points, especially well seen for m = 4.. 6. Larger m 
corresponds to smaller angular size of the sound circles. 
This happens because we have not fitted the freezeout 
temperature and time Tf to these data: decreasing the 
former and increasing the latter one can certainly get 
better fit. We have not done so because in any case our 
calculation is done for conformal matter with fixed speed 
of sound and e/T^, and cannot accurately describe the 
real collisions anyway. 



E. The location of the perturbation 

So far we have demonstrated some qualitative features 
of the one-body spectrum and two-body correlations re- 
sulting from a local perturbation, selecting one typical 
location. In this section we provide further detail on the 
modifications of the Green function we calculated on the 
location of the initial hot spot. Since we only consider 
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FIG. 10: (Color online) Top: The two-pion distribution in 
arbitrary units as a function of azimuthal angle difference Aij) 
(rad), for r =2(blue large dash),3(brown dash-dot), 4.1 (red 
solid line) fm. Bottom: The two-pion distribution in arbitrary 
units as a function of azimuthal angle difference A(j> (rad), 
for r =4.1(the same red solid line), 4. 7 (green small dash), 5. 5 
(black dash-dot-dot) fm. All plots are for the same value of 
the viscosity-to-entropy ratio rj/s — 0.134 

central collisions, by "location" we mean the radial posi- 
tion of the "hot spot" . As shown in Fig 10 changing the 



location of the spot visibly affects the quantitative shape 
of the two-particle correlation as well as the power spec- 
trum FigjlT] When the spot is located near the center 
of the fireball, the two particle correlation presents only 
one peak located at A0 = 0, and no structure on the 
away side. The characteristic two peaks appear when 
the initial perturbation is located not too close to the 
center (r ~ 3 — 5 fm). 

Furthermore, as one can see, the amplitude of the mod- 
ulation decreases in this case. This happens not because 
of a change of the hot spot amplitude (which is the same 
in all cases), but because of the (partial) cancellation be- 
tween hydro perturbations for velocities of the first type 
(in the sound wave) and the second type (extra radial 
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FIG. 11: (Color online) The(red) solid, (blue) dashed, (green) 
dash-dotted and (black) dotted curves correspond spectral 
distributions obtained for initial perturbations located at 
r = 3,4.1,4.7 and 5.5 fm, respectively, for 77/s = 0.134. 

flow stemming from the modification of the freezeout sur- 
face). As we have discovered, the very sign of the projec- 
tion of the former on the radial direction depends on the 
initial position of the perturbation. For perturbations 
located near the center of the fireball it is positive, but 
as the "hot spot" gets located at larger r, it decreases 
becoming negative till it gets as large as the second one 
and cancels it, when the "hot spot" is located at the very 
edge of the fireball. 

In Fig {11] it is possible to see how the change in the ra- 
dial position of the initial perturbation affects the power 
spectrum. Its general features remain unaltered, present- 
ing maxima and minima in all cases, which decrease for 
larger values of m due to viscosity. The figure shows that 
there is some shift with r in the position of the maxima 
and minima. 

In order to compare our results with the experimen- 
tal data, it would be necessary to average over different 
initial perturbations, using probability distributions for 
their locations and amplitudes. Since the minima for the 
different locations do not precisely match, in an aver- 
aged case a minimum would still be present, but it would 
not be as pronounced as in the case of an individual ini- 
tial perturbation, the whole shape of the power spectrum 
would be smoother, with no sharp dips. In principle, very 
precise data can potentially be used to infer some infor- 
mation about the perturbation distribution in r. Such 



averaging is deferred to the subsequent works, as it would 
require a particular model for the initial state. It can be 
the Glauber model (we discussed in our previous paper) 
or some models including the saturation phenomenon. 



V. SUMMARY AND FINAL COMMENT 

By calling this work "the second act of hydrodynam- 
ics" we emphasize the huge progress made in the field. 
From measuring the mean velocity of matter and the 
mean ellipticity a decade ago, the first evidences for col- 
lective fiow, we now have data providing up to the 9- 
th harmonics of it. With many theory results, some of 
them in this work, we also now have an understanding 
of how perturbations behave as m grows. In short, the 
answer is that they are acoustic oscillations, with certain 
m-dependent oscillation frequencies and dampings. We 
have found that, like in the Big Bang, rotating phases 
at the freezeout generate minima and maxima. Remark- 
ably, experimental data provide the first hints for the 
minimum and the second maximum. 

The rather intricate shape of the two-particle correla- 
tions as a function of A0 is very similar to the results of 
our calculation of the Green function from a local source. 
But we would like to mention, as a parting comment, that 
the questions: Do the sound circles exist in reality, or is 
it just a mathematical tool ? Are different harmonics co- 
herent or not? are still unanswered an they represent the 
next challenge for the field. A way to figure this out is 
explained in our previous paper [T]: one should measure 
the t/iree-particle correlation functions, and look for the 
"resonances" between 3 harmonics related by the "tri- 
angular" condition mi -|- 1712 + TO3 = 0, or by the two- 
particle correlations with respect to reaction plane (for 
non-central collisions). 
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